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KELSEA 

AN  INTERACTIVE  COMPUTER  CODE 
FOR 

KELVIN  AND  RANDOM  AMBIENT  SEA  WAVES 


I.  INTRODUCTION 

In  many  ocean  engineering  applications  it  is  of  interest  to  generate  a  spatial  realization  of  a  given 
random  sea  state  from  an  experimentally  derived  or  analytical  wave  spectrum.  Such  two  dimensional 
spatial  realizations  can  be  used  to  model  the  wave  elevations  around  a  large  offshore  structure  or  to 
model  the  interaction  of  these  waves  with  ocean  currents  and/or  ship  waves. 

There  is  a  reasonably  large  number  of  studies  in  the  literature  on  the  general  approach  for  gen¬ 
erating  spatial  realizations  of  random  ambient  sea  waves  as  a  double  sum  of  sinusoidal  components 
from  a  given  energy  spectrum.  The  general  approach  is  relatively  simple  and  it  is  briefly  discussed  in 
the  following  section  of  this  paper.  Detailed  descriptions  of  the  general  approach  can  be  found  in  [1]- 
[3].  Although  the  basic  approach  has  been  well  documented,  its  manner  of  application  has  varied 
widely.  For  example,  the  number  and  assigned  location  of  the  components  required  to  obtain  a  random 
ambient  sea  state  are  not  always  clearly  stated  and  a  trial  and  error  process  is  often  used  to  obtain  the 
wave  elevation  signal. 

The  present  paper  addresses  the  specific  ways  of  implementing  the  generalized  approach  and 
theory  for  generating  a  random  ambient  sea,  based  on  the  approach  given  in  [3j.  It  has  been  shown  in 
this  reference  that  relatively  few  components  are  required  to  obtain  sea  states  which  visually  and  statist¬ 
ically  approximate  those  for  a  random  sea. 

The  Kelvin  waves  generated  by  a  ship  traveling  at  constant  forward  speed  are  obtained  by  using 
the  complex  elevation  spectra,  which  can  be  calculated  by  the  computer  program  KELVIN  14] .  The 
particular  method  used  to  calculate  the  Kelvin  wave  elevations  over  a  rectangular  grid  is  described  in 
some  detail  and  the  linear  superposition  of  the  Kelvin  and  ocean  waves  is  also  discussed. 

The  mathematical  model  has  been  implemented  in  a  computer  program  which  is  described  in  the 
last  section.  Examples  of  selected  cases  of  random  ambient  sea  states  are  given  as  representative  exam¬ 
ples  in  the  same  section.  Furthermore,  examples  are  also  given  for  the  superposition  of  Kelvin  waves 
with  ambient  sea  waves. 

II.  MATHEMATICAL  MODEL 

The  development  of  random  sea  state  realizations  is  based  on  expressions  which  involve  a  triple 
sum  over  the  wave  number  k,  frequency  w,  and  direction  9.  If  one  neglects  the  nonlinear  interactions 
between  wave  components,  the  following  expression  for  the  wave  elevation  tj(x,  z)  can  be  derived 
from  the  disperson  relation  between  wave  number  and  frequency 

N  M 

t)(x,  z)  =*  II  AtJ  cosU,  Or  cos0,  +  z  sin0y)  +  <£,,]  (1) 

.-i /-I 

where  x  and  z  are  respectively  the  coordinates  parallel  and  normal  to  the  predominant  wave  direction, 
A,,  is  the  amplitude  of  the  (i,  j) th  wave  component,  S(u»,,  9,)  is  the  value  of  the  energy  spectrum,  9 
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is  the  wave  direction  measured  from  the  X-axis,  k  =  u>}/g  is  the  wave  number  for  deep  water  and  0„ 
is  a  random  phase  angle  uniformly  distributed  between  9  and  2n.  For  the  case  where  S(w,  0)  is  the 
energy  spectrum,  Ai}  is  given  by 

Au  =  ^2SU^~0/TK^9  (2) 

where  Aw  and  A 9  are  increments  of  w  and  9 ,  respectively.  This  term  should  be  treated  correctly  in 
order  to  avoid  errors  of  a  factor  of  two  in  the  wave  elevations.  In  the  present  study  all  spectra  are 
given  in  terms  of  energy,  which  is  consistent  with  Eqs.  (1)  and  (2)  and  it  is  assumed  to  have  the  com¬ 
monly  used  separable  form 

S(w,  9)  =  SuU)  G(0)  (3) 

where  S<„(w)  is  the  frequency  spectrum  and  G(0)  is  the  spreading  function  which  must  satisfy  the 
equation 


G(9)d0  =  1 

-it 


Three  widely  used  frequency  spectra  are  presently  considered:  the  Pierson-Moskowitz,  the 
Jonswap  and  the  Neumann.  A  fourth  one,  the  Bretschneider,  which  is  a  subset  of  the  Pierson- 
Moskowitz,  is  also  included  in  the  computer  program  as  a  special  case.  In  their  most  general  form  the 
first  two  are  functions  of  the  wind  velocity  and  they  are  comparable  to  the  Neumann  spectrum  if  only 
forms  appropriate  to  fully  developed  seas  are  considered.  To  conform  with  this  requirement  an  expres¬ 
sion  for  the  three  spectra  may  be  written  in  the  following  form 

Sm  =  e~[BI •  3.3’  (5) 

Qjm 

where  3.3  is  a  mean  enhancement  factor  and  the  values  of  the  other  parameters  for  each  of  the  spectra 
are  given  in  Table  1,  where  w0  is  the  peak  frequency  given  by 


w0  =  .8771632  (*/U) 

and  cr  is  equal  to  0.07  for  w  <  w0  and  equal  to  0.09  for  w  >  w0. 

Table  1 


Pierson- 

Moskowitz 

Jonswap 

Neumann 

A 

0.008 1*G**2 

0.0081*G**2 

33.1V/4 

B 

0.74*G**4/U**4 

0.74*G**4/U**4 

2*G**2/U**2 

m 

5 

5 

6 

n 

4 

4 

2 

q 

0 

[— <0— Wq)2/  (2o-2t«J^) 

0 

A  number  of  similarities  exist  among  the  three  frequency  spectra  as  well  as  some  distinct  differ¬ 
ences.  The  Pierson-Moskowitz  and  Jonswap  spectra  are  similar  but  with  a  significant  exception  that  the 
Jonswap  is  more  sharply  peaked  in  the  region  of  the  peak  frequency  w.  The  Neumann  spectrum  differs 
from  the  other  two  in  that  the  values  for  the  parameters  m  and  n  are  different.  The  Bretschneider 
spectrum  has  the  same  values  as  the  Pierson-Moskowitz  for  the  m,  n  and  <7,  but  it  differs  in  the  values 
for  A  and  B. 

Several  types  of  expressions  can  be  found  in  the  literature  for  the  spreading  function  G(9 V  A 
relatively  simple  yet  reasonably  accurate  expression  for  this  function  is  the  widely  used  cosine  squared 
function  given  by 
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where  v  is  the  coefficient  which  determines  the  range  of  0  over  which  the  energy  of  the  wind  is  spread. 
The  simplest  and  most  popular  form  for  (7(0)  is  to  let  p  —  1.  In  this  form,  the  wind  is  assumed  to 
generate  waves  in  all  directions  which  are  within  ±90  degrees  from  the  wind  direction.  Values  of  v  — 
3/2  and  v  —  3/4  correspond  to  spreadings  of  ±60  degrees  and  ±120  degrees  respectively. 

The  definitions  and  specification  of  the  frequency  spectra  and  spreading  function  are  required  for 
the  calculation  of  the  wave  amplitude  coefficients  A0 ,  which  are  given  by  Eq.  (2).  The  last  parameters 
requi.:d  by  Eq.  (1)  are  the  values  of  the  x  and  z  coordinates.  These  are  specified  by  defining  a  rec¬ 
tangular  grid  over  the  sea  surface,  which  is  described  in  the  following  section.  Then,  at  each  of  the  grid 
points  the  sea  surface  elevations  are  calculated  by  using  Eqs.  0),  (2),  (S)  and  (7).  These  equations 
together  with  the  parameters  of  Table  1  have  been  implemented  in  the  computer  program  KELSEA  for 
calculating  the  wave  elevations  for  random  ambient  sea  waves. 

III.  SUPERPOSITION  OF  AMBIENT  AND  KELVIN  WAVES 

The  Kelvin  wave  elevations  f  (x,  z)  are  obtained  by  taking  the  inverse  Fourier  transform  of  the 
complex  wave  amplitude  spectrum  A(u)  calculated  by  the  program  KELVIN  (4). 

*  (x,  z)  =  4?  T-Re  f~  A  iu)e-,k^m)du  (8) 

where  u  =  k:/k0  is  the  dimensionless  lateral  wavenumber,  s  =  kx/k0  is  the  dimensionless  longitudinal 
wavenumber,  k0  =  g/  U** 2  is  the  fundamental  wavenumber  and  Re  represents  the  real  part  of  the  com¬ 
plex  integral.  The  Kelvin  waves  are  calculated  for  a  coordinate  system  fixed  to  the  ship  and  with  its 
origin  at  the  forward  perpendicular.  The  x  and  z  coordinates  are  in  the  longitudinal  and  lateral  direc¬ 
tion  respectively  and  on  the  horizontal  plane.  The  integral  in  Eq.  (8)  is  calculated  for  evenly  spaced 
values  of  the  dimensionless  wavenumber  u.  The  values  of  x  and  z  are  chosen  to  lie  in  the  rectangular 
grid,  with  equal  spacing  in  each  direction.  The  extent  of  the  grid  in  the  lateral  z  direction  is  equal  to 
2X0,  (±  X0)  on  either  side  of  the  x  axis,  where  X0  is  the  fundamental  wavelength  given  by 

2 it  _  4tt  V 2 
k<>  g 


Xn  = 


(9) 


Noting  that  the  cusp  line  of  the  Kelvin  waves  form  an  angle  of  19.5  degrees  with  the  x  axis, 
required  extent  of  the  grid  in  the  x  direction,  XM,  is  simply  given  by 


the 


X*  = 


tan(I9.5) 


(10) 


The  computer  program  KELSEA  has  the  options  of  generating  the  ocean  and  Kelvin  wave  eleva¬ 
tions  separately  or  a  superposition  of  the  two.  In  performing  the  superposition,  it  is  of  interest  to  note 
that  the  ocean  waves  resulting  from  Eq.  (1)  are  generated  for  a  coordinate  system  which  is  fixed  in 
space  while  the  Kelvin  waves  are  generated  for  a  coordinate  system  which  is  fixed  to  the  ship.  In  the 
general  case  where  the  wave  elevations  are  functions  of  time,  they  would  differ  in  the  two  coordinate 
systems  since  the  wave  frequencies  encountered  by  a  moving  ship  are  different  from  those  experienced 
by  a  stationary  one.  However,  in  the  present  case,  time  is  not  a  variable  since  the  Kelvin  waves 
expressed  by  Eq.  (8)  are  steady  in  the  moving  coordinate  system  and  the  ocean  waves  expressed  by  Eq. 
(1)  represent  an  image  of  the  ocean  surface  at  a  particular  instant  of  time,  say  t  =  0.  Consequently, 
direct  superposition  of  the  waves  is  permissible,  without  the  need  to  account  for  the  difference  in  the 
velocities  of  the  two  coordinate  systems. 

Finally,  the  case  of  ambient  waves  with  the  wind  direction  at  an  angle  0„  with  the  x-axis,  is  con¬ 
veniently  treated  by  simply  substituting  0'  -  0  +  0U  in  Eq.  (1).  Example:,  of  ambient  ocean  waves  for 
different  wind  directions  are  presented  in  the  next  section. 
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IV.  COMPUTER  PROGRAM 


The  mathematical  model  discussed  in  the  previous  sections  has  been  incorporated  into  the  com¬ 
puter  program  KELSEA.  The  major  tasks  performed  by  this  program  are  discussed  in  this  section  as 
well  as  the  structure  of  the  KELSEA  program.  Furthermore,  the  interactive  capabilities  of  the  program 
are  demonstrated  by  a  typical  example  run  and  calculated  results  for  certain  input  parameters  are  also 
given. 

The  KELSEA  program  is  designed  to  perform  three  functions:  a)  It  performs  the  random 
ambient  sea  state  calculations  according  to  specifier  iut  parameters,  b)  It  computes  the  free  surface 
elevations  for  Kelvin  waves  from  a  pre-calculated  >r  .tral  file,  and  c)  It  superimposes  the  surface  eleva¬ 
tions  for  a  Kelvin  wave  and  the  ambient  sea  wav  ,.  The  superposition  is  a  linear  point  to  point  super¬ 
position  and  should  not  be  interpreted  as  an  interaction  between  the  sea  waves  and  the  Kelvin  wave. 
The  code  is  written  in  FORTRAN-77  and  uses  menus  and  inquiries  to  enhance  its  use  and  it  is 
currently  operating  in  a  number  of  computer  systems  such  as  the  HP-1000,  HP-9000  and  VAX  11/780 
systems.  The  modular  structure  of  the  computer  program  is  shown  in  Fig.  1. 


Fig.  1  —  KELSEA  computer  program  diagram 


KELSEA  is  programmed  to  perform  three  independent  functions.  Each  function  is  directed  by  a 
program  menu.  The  first  menu  item  is  the  Kelvin  wave  heights  calculations.  Illustrated  below  is  the 
program  flow  diagram  for  calculating  the  free  surface  elevations  of  the  Kelvin  waves. 
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This  task  requires  the  input  of  a  Kelvin  spectral  file,  which  is  created  from  the  KELVIN  program 
First,  the  free  surface  grid  is  generated  in  the  SURFACEGRID  routine.  Then  the  surface  elevations 
for  the  Kelvin  waves  are  calculated  in  the  KELVINWAVE  routine  and  the  resulting  output  data  are 
stored  in  a  file  for  plotting. 

The  next  task  on  the  menu  is  the  sea  state  calculations.  Its  input  data  is  program-generated  based 
on  user  specified  parameters  and  spectrum  type.  It  is  important  to  note  here  that  for  every  Kelvin 
spectral  file  there  must  be  a  sea  state  calculation  corresponding  to  the  same  surface  grid  in  order  to 
correctly  superimpose  the  calculations.  The  reason  is  that  velocity  will  vary  for  every  Kelvin  spectral 
file  and  the  grids  will  be  dimensioned  differently.  Once  the  surface  grid  is  created  and  the  sea  state 
elevations  are  computed  in  the  COSQSP  routine,  the  results  are  written  into  a  file  by  the  WRITE  FILE 
routine.  The  following  diagram  shows  the  program  flow  of  the  sea  waves  calculations. 
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The  third  task  performed  by  the  KELSEA  program  is  the  linear  superposition  of  the  Kelvin  wave 
and  sea  state  calculations. 
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This  task  requests  for  two  input  files:  the  Kelvin  wave  surface  elevations  data  file  and  the  sea  wave 
elevations  data  file.  These  two  files  must  be  created  before  the  linear  superposition  of  the  two  can  be 
performed.  The  previous  tasks  explained  above  generate  the  output  files  that  can  V  used  as  the  input 
here.  The  computations  take  place  in  the  KELVIN  AN D_SE A  routine  and  the  results  are  written  into 
a  file  by  the  WRITE  FILE  routine.  Any  of  these  output  files  can  be  plotted  by  the  PLOT3D  program 
(5]  as  three  dimensional  surfaces. 

The  three  tasks  for  the  KELSEA  program  are  controlled  by  the  program  main  menu. 

The  user  is  requested  to  make  a  selection  from  the  program  menu  as  illustrated  below. 
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IV. 1  Kelvin  Wave  Calculations 


To  run  the  Kelvin  wave  calculations  enter  the  corresponding  number. 

**  Specify  Menu  Item  :  1 


The  program  will  inquire  for  the  name  of  the  Kelvin  spectral  datajjle't6"be  read  and  stored  for 
calculations. 


******************«*****************^4«********************* 
****  c ,u„  n  . . £;i*,  f-„„.  «*»* 

»**»  .  u ; d  .. .. . - ,.,;n  d  ....  a  .  u i ...... .  ra . .,  »*** 

**♦* 


Specify  the  Name  of  the  Data  File  from  which 
this  Program  will  Read  the  Input  Data 
File  Name  :KELVIN.DAT 


The  description  of  the  Tile  contents  are  printed  on  the  screen  as  follows 


Data  Read  from  File 

:  KELVIN.DAT 

Velocity  (Ft. /Sec.)  = 

30.000 

Characteristic  Wave  #  = 

.04 

Y  Dist  to  P  Symmetry 

0.00 

Displace  men  t  (Tons)  = 

34461.18 

Waterline  Length  (Ft.)  = 

300.00 

Wet  Surf  (Ft.**2) 

1041.00 

L  Center  Buoy.  (Ft.)  = 

-136.56 

#  of  Pos.  Wave  Nos.  = 

350 

2  Hulls?  (Y=1/N  =  0) 

0 

The  next  input  required  by  the  program  is  the  number  of  points  for  the  rectangular  grid  in  the  X 
and  Z  directions  and  the  starting  X-coordinate  for  grid.  Currently  the  axes  origin  is  at  the  forward  per¬ 
pendicular  with  the  X-axis  along  the  ship  length  taken  as  negative  and  with  the  starting  position  equal 
to  a  ship  length.  The  user  has  the  option  of  using  the  default  values  given  or  define  new  values  for  the 
number  of  points.  In  this  case,  enter  the  default  values. 

*************************************************************** 

****  Specify  the  Number  of  Points  for  the  Surface  **** 

****  Grid  in  both  X-  and  Z-  directions. 

****  Default  values  are  NXPOINT  =  61,  NZPOINT  =  41 
****  NX-POINT  =61 
****  NZ-POINT  =41 


****  Specify  the  starting  value  in  the  X  direction 
****  for  the  rectangular  grid  in  (ft).  The  default 
****  value  is  equal  to  the  Ship  Length.  XMIN  = 


Once  the  surface  grid  has  been  defined,  the  program  computes  the  free  surface  heights  for  the 
Kelvin  wave.  After  the  calculation  process  is  completed,  the  next  step  is  to  write  the  results  to  a  file. 
The  program  indicates  the  type  of  results  which  are  to  be  saved  as  follows 

>>>>>>>>>>  KELVIN  WAVE  OUTPUT  FILE  <<<<<<<<<< 

and  it  requests  for  a  name  of  the  output  file  to  save  the  results. 


.  ^  _ 


«■  * 


*********************************  *************************** 

****  Specify  the  Name  of  th.  Data  File  For  Storing 
****  the  Output  Results  for  Plotting 
****  File  Name  :KELPLOT.DAT 


The  resulting  information  about  the  file  just  created  is  written  on  the  screen,  which  concludes  the 
task  processing  for  the  Kelvin  wave  calculations,  and  the  program  returns  to  the  main  menu. 

Data  Written  to  File  :  KELPLOT.DAT 

NXPOINT  =  61 
NZPOINT  =  41 


IV. 2  Sea  State  Calculations 

The  program  menu  appears  on  the  screen  for  the  user  to  make  another  selection.  To  calculate  the 
sea  state  elevations  enter  the  corresponding  number  (2)  to  the  menu  inquiry. 

The  following  notice  is  displayed  on  the  screen  reminding  the  user  that  a  corresponding  Kelvin 
file  must  exist  for  the  same  grid  dimensions  and  velocity. 

************************************************* 

....  NOTICE:  »•'* 

****  This  task  requires  that  a  Kelvin  *"* 

....  spectral  file  has  been  read  to  **** 

****  determine  the  parameters  needed  **** 

****  for  the  surface  grid  calculations.  ’*'* 

****  If  a  file  is  not  available  you  will  **•* 

****  be  asked  to  specify  the  dimensions  ’*** 

****  of  the  physical  space  for  the  **** 

****  sea  slate  calculations.  **** 


The  next  request  by  the  program  is  the  specification  of  the  spectrum  type  One  of  the  options 
given  below  is  selected  by  entering  the  corresponding  number. 

**  specify  the  Spectrum  Type 

**  for  the  Sea  Slate  Calculations 

**  Input  (I)  for  P1ERSON-MOSKOWITZ 

**  (2)  forJONSWAP 

**  (3)  for  NEUMANN 

*•  (4)  lor  BRFTSC1INEIDER 

**  Specify  the  Approiate  Spectrum  .  1 

The  inquiries  below  are  prompted  one  at  a  time  in  order  to  specify  the  appropriate  value. 

**  For  60  Deg  ,  Spreading  index  =  1.50** 

**  For  90  Deg.,  Spreading  index  =  1.00** 

**  for  120  Deg.,  Spreading  index  =  0  75  ** 

**  Specify  the  Spreading  Index  (ge  51  :  0.5 

**  Specify  the  Seed  Number 

**  for  Random  Phase  Angie  < ISFED  fiE  1 1  :  1 


>.V„ 


Specify  Wind  Speed  [S.S.  3  =  16  Kts]  :  16 


**  Specify  the  Wind  Direction 
**  with  respect  to  the  X-axis  in  deg.  :  45 

**  Specify  the  Number  of  Wave 
**  Frequencies  (1  to  101)  :  10 

**  Specify  ZERO  (0)  for  Equal  ASW 
**  or  ONE  (1)  for  Equal  Delta  Freq.  :  0 

The  next  input  required  by  the  program  is  the  number  of  points  for  the  surface  grid  calculations 
and  the  starting  value  in  the  X  direction. 


****  Specify  the  Number  of  Points  for  the  Surface  **** 

****  Grid  in  both  X-  and  Z-directions.  **** 

****  Default  Values  are  NX_POINT  =  61,  NZ_POINT  =  41  ***• 

****  NX-POINT  =61 
***•  NZ-POINT  =41 


****  Specify  the  starting  value  in  the  X  direction  **** 

****  for  the  rectangular  grid  in  (ft).  The  default  **** 

****  value  is  equal  to  the  Ship  Length.  XM1N  = 

If  the  user  were  to  calculate  the  sea  wave  elevations  only,  then  the  inquiries  below  would  be 
made.  Since  in  the  current  example  a  KELVIN  file  was  read,  the  program  will  use  that  information  for 
the  actual  dimensions  of  the  surface  grid.  Remember  that  the  surface  grid  for  the  Kelvin  and  sea  wave 
elevations  must  be  of  the  same  dimension  and  relating  to  the  same  ship  velocity  in  order  to  make  a 
rational  superposition  of  the  resulting  calculations. 

**•♦*♦*******•******•*******«*•••*«**•****•*•*•••«**•**•*« 

****  Specify  the  Value  for  the  Rectangular  Grid  **** 

****  in  the  X  direction  in  ft.  (  XMAX  )  **** 

****  XMAX  = 


****  Specify  the  Value  for  the  Rectangular  Grid  **** 

****  in  the  Z  direction  in  ft.  (  -ZMIN  ,  +  ZMAX  ) 

****  ZM AX  = 

The  above  two  inquiries  can  be  used  if  one  wishes  to  generate  an  ambient  sea  state  independent 
from  a  Kelvin  wave. 

The  free  surface  grid  and  elevations  are  calculated  for  the  sea  state  task.  The  output  message 
below  indicates  to  the  user  which  results  have  been  computed  and  to  be  written  to  a  file. 

>>>>>>>>  AMBIENT  SEA  WAVE  OUTPUT  FILE  <<<<<<<< 

The  user  is  to  input  the  file  name  for  the  sea  wave  results. 

****  Specify  the  Name  of  the  Data  f  ile  For  Storing  **** 

****  The  Output  Results  for  Plotting  **** 

**•*  File  Name  :SEAPLOT.DAT 


IV. 3  Superposition  of  the  Kelvin  and  Sea  Wave  Elevations 


The  program  returns  control  to  the  main  Menu  for  a  selection  to  be  made.  Enter  the  correspond¬ 
ing  number  (3)  for  the  next  task. 


There  will  be  a  message  printed  on  the  screen  reminding  the  user  that  before  the  program  contin¬ 
ues,  both  data  files  for  the  Kelvin  and  sea  wave  elevations  must  exist,  for  the  linear  superposition  to  be 
processed.  Answer  the  inquiry  and  continue. 


**** 
**** 
**** 
**♦* 
**** 
**** 
**♦* 
*  *  *  * 
**** 
**** 
**** 
*  *  *  * 


NOTICE  : 

This  task  requires  that  the  Kelvin 
wave  and  the  sea  wave  elevation 
data  files  have  been  created. 


Enter  (Y)  if  the  files  exist  and 
continue  with  the  calculations. 

Or,  enter  (N)  if  the  files  have  not 
been  created  and  return  to  the 
program  menu  to  create  the  files. 


♦  *** 
**** 
*  *  *  * 
**** 
**** 
**** 
**** 
**** 
**** 
#*** 
**** 


ENTER  Y/N  :Y 


The  program  will  inquire  for  the  names  of  these  two  data  files  and  it  will  check  for  consistency  of 
the  two  surface  grids.  If  the  surface  grids  corresponding  to  the  two  specified  file  are  not  the  same  an 
error  message  will  be  printed  on  the  screen. 

******************************************************** 

****  Specify  the  Name  of  the  Input  File  for  the  **** 

***♦  KELVIN  WAVE  Results 
*♦**  File  Name  KELPLOT.DAT 


**** 
**♦* 


****  Specify  the  Name  of  the  Input  File  for  the 

****  SEA  STATE  Results 

****  File  Name  :SEAPLOT.DAT 


The  output  message  below  indicates  to  the  user  which  output  results  are  to  be  stored  in  a  file  for 
later  use. 


>>>>>>>>  KELVIN/SEA  WAVE  OUTPUT  FILE  <<<<<<<< 


**** 
**** 


****  Specify  the  Name  of  the  Data  File  For  Storing 
****  The  Output  Results  for  Plotting 
***'  “  •••  Mtme  :KELSEA.DAT 


Once  again  the  program  menu  appears  on  the  screen.  All  the  tasks  have  been  processed  and  the 
user  can  exit  the  program  by  entering  the  corresponding  number  (4)  to  the  inquiry. 


V.  NUMERICAL  RESULTS 


The  three  output  files  created  by  the  KELSEA  program  can  be  plotted  as  three  dimensional  sur¬ 
faces.  For  the  purpose  of  providing  some  examples  of  numerical  results  in  this  paper,  a  spectral  input 
file  will  be  used.  This  lilc  has  been  generated  by  the  program  KELVIN  (4) ,  which  calculates  the  Kelvin 


v.v 


.  -  .  A< 


.V.v 
V.V,V 
VVX' 


V.  V 


» 

:v 

-V.V. 


•yv. 


•y- 


•r.  1 
V.V*' 
*  ,  . 

/vV 

k:' 


•r.  V.  •*. 

'•“.v.v-,-. 

y.v. 


,>Vv 

s  • " 


J  ,■  s*  *  ,■ 


*  «T.  *  f  •  , 


V- 


wave  for  a  specific  ship  hull  configuration.  In  this  particular  example  the  spectal  file  corresponds  to  a 
mono-hull  ship  configuration.  The  relevant  parameters  from  the  input  file  are  given  in  Table  2. 

The  resultinS  output  files  for  the  Kelvin  wave  elevations  are  designated  as  KELWAV1.DAT  and 
KELWAV2.DAT  respectively.  Similarly  the  output  files  for  the  sea  wave  elevations  are  designated  as 
SEAWAV11.DAT,  SEAWAV12.DAT,  SEAWAV21.DAT  and  SEAWAV22.DAT.  The  input  parame¬ 
ters  used  for  generating  these  four  files  are  given  in  Table  3. 


Table  2 


Input  File  :  KELVIN1.DAT 

Velocity  (Ft./Sec.)  =  31.45 

Characteristic  Wave  #  =  .03 

Y  Dist  to  P  Symmetry  =  0.00 

Displacement  (Tons)  =  2599.60 

Strut  Length  (Ft.  )  =  300.00 

S  Wet  Surf  (Ft.**2)  =  15274.50 

L  Center  Buoy.  (Ft.)  =  -157.825 

#  of  Pos.  Wave  Nos.  —  300 

2  Hulls?  (Y=1/N=0)  0 


Input  File  :  KELVIN2.DAT 

Velocity  (Ft./Sec.)  =  31.45 

Characteristic  Wave  #  =  .03 

Y  Dist  to  P  Symmetry  =  80.00 

Displacement  (Tons)  =  5199.20 

Strut  Length  (Ft.  )  =  300.00 

S  Wet  Surf  (Ft.**2)  =  30549.00 

L  Center  Buoy.  (Ft.)  =  -157.825 

#  of  Pos.  Wave  Nos.  =  300 

2  Hulls?  (Y  =  1/N=0)  1 


Table  3 


Input  file 
Output  Files 


KELVIN1.DAT 

KELWAV1.DAT 

SEAWAV11.DAT 

SEAWAV12.DAT 


Input  file 
Output  Files 


KELVIN2.DAT 

KELWAV2.DAT 

SEAWAV21.DAT 

SEAWAV22.DAT 


Sea  State  case  :  ONE 

Spectrum  Type 
Spreading  Index 
Seed  Number 
Wind  Velocity 
Wind  Direction 
Number  of  wave  Freq. 


Sea  State  case  :  TWO 

Spectrum  Type 
Spreading  Index 
Seed  Number 
Wind  Velocity 
Wind  Direction 
Number  of  wave  Freq. 


The  files  SEAWAV11.DAT  and  SEAWAV21.DAT  correspond  to  case  one  for  the  sea  state. 
Similarly  SEAWAV12.DAT  and  SEAWAV22.DAT  correspond  to  case  2  for  the  sea  state.  The  wind 
velocity  of  16  knots  corresponds  to  a  Sea  State  of  (3). 

The  numerical  results  from  the  six  output  files  are  presented  in  Figs.  2-7  as  follows 


KELWA1.DAT, 

Fig.  2 

KELWAV2.DAT, 

Fig.  5 

SEAWAV11.DAT, 

Fig.  3 

SEAWAV21.DAT, 

Fig.  6 

SEAWAV12.DAT, 

Fig.  4 

SEAWAV22.DAT, 

Fig.  7 
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The  superposition  of  each  of  the  Kelvin  wave  elevation  files  with  each  of  the  four  sea  state  files 
will  result  in  the  following  output  files: 

KELWAV1.DAT  +  SEAWAV11.DAT  —  KELSEA111.DAT,  Fig.  8 

KELWAV1.DAT  +  SEAWAV12.DAT  —  KELSEA112.DAT,  Fig.  9 

KELWAV1.DAT  +  SEAWAV21  DAT  —  KELSEA121.DAT,  Fig.  10 

KELWAV1.DAT  +  SEAWAV22.DAT  —  KELSEA122.DAT,  Fig.  11 

KELWAV2.DAT  +  SEAWAV11.DAT  —  KELSEA211.DAT,  Fig.  12 

KELWAV2.DAT  +  SEAWAV12.DAT  —  KELSEA212.DAT,  Fig.  13 

KELWAV2.DAT  +  SEAWAV21.DAT  —  KELSEA221.DAT,  Fig.  14 

KELWAV2.DAT  +  SEAWAV22.DAT  —  KELSEA222.DAT,  Fig.  15 
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The  plotted  results  of  these  eight  files  are  presented  in  Figs.  8  through  15.  These  results  have 
been  plotted  with  the  PLOT3D  computer  program. 

VI.  CONCLUSIONS 

The  mathematical  model  for  calculating  ambient  sea  waves  has  been  presented  in  this  report.  The 
modeling  and  calculations  of  the  Kelvin  wave  elevations  and  their  superposition  onto  the  ambient  sea 
waves  also  has  been  presented. 

The  computer  program  implementing  these  concepts  has  been  described  and  the  program  struc¬ 
ture  has  been  discussed.  A  number  of  numerical  examples  have  been  given  in  order  to  illustrate  the 
program  capabilities. 

The  computer  program  has  been  developed  in  such  a  way  that  it  is  simple  to  operate  and  the  user 
does  not  need  to  refer  to  manuals  for  exercising  the  program.  Furthermore,  the  input  parameters  for 
the  sea  wave  calculations  add  flexibility  to  the  types  of  Sea  States  that  one  can  generate 
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Figure  15 
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Appendix 

KELSEA  Program  Listing 


PROGRAM  KELSEA 


* 

*  Revi sion  List: 


*  --Date--  --by--  --Description  — 

* 

*  Jun.  1985  -  First  Revised  Version  for  H-P  Computers,  Ver .  1.0 

*  Sep.  1985  -  Final  Version  for  H-P  Computers,  Ver.  1.2 

*  Dec.  1985  -  VAX  Version,  Ver.  1.2V 

* 


* 

*  *  * 
*****  Free  Surface  Heights  Calculation  for  a  Kelvin  Wave 
*****  and  Sea  State  Ca 1 cu 1  a t ions  . 

★  * 

COMMON/BLOCK1 /  X_POINT(  101  )  , Z_POINT(  101  )  ,Y_PO INTI  101  ,101  ) 

COMMON/ BLOCK 2/  CGG( 100  0 )  ,CGM( 10  00  ) 

COMMON /BLOCK  3/  U_KELVIN(  1000  )  ,V_KELV1N( 1  000 )  , W_KELVIN(  1000 ) 
COMMON/ BLOCK4/  NX_POINT , NZ_POINT , NY_POINT , DX , DZ , GR 1 D_S I ZE 
COMMON/ BLOCKS /  VEE , AK , YSYM , DISPT , XST , ZST , AREAS , CB , N_POI NTS , NSYM 
COMMON/ BLOCK 6 /  WM(101),AZD(101),WZD(101)  ,  ENF ( 101  )  , ENl 161 ) 

COMMON/ BLOCK7/  NSP , FNXP , I  SEED , WVKT , THSTD , TZ , NF , 1 EDW 
DIMENSION  LBUF(IOOO) 

CHARACTER*  20  FILE1 ,FILE2 
COMPLEX  CGG , CGM 


*  * 

CALL  LGBUF( LBUF  ,  1000  ) 

** 

GR  I  D_S  I  ZE  --  0.0 

?  Default 

value 

NX_POI NT  =  61 

!  Default 

va  1  ue 

NZ_PO I  NT  =  41 

!  Default 

value 

*  * 


10  CALL  MENU ( I DMENU ) 

*  * 

GO  TO  (100,200,300,500)  I DMENU 
100  CONTINUE 

*  * 

*****  Read  Spectral  File  for  Height  Calculations. 

*  * 

CALL  READ_FILE( FILE1  ) 

*  * 

WRITE( 1 , ' ( / ,5X, “Data  Read  from  File  :  ",A20))  F1LE1 

WRITE (1, ' (/,5X," Velocity  (Ft. /Sec.)  =",F10.2,/ 

+  5X , “Cha r ac t er l s t l c  Wave  #*",F10.2,/ 

+  5X,"Y  Dist  to  P  Symmetry  =",F10.2  )')  VEE, AK, YSYM 

WRITE ( 1 ,'(/, 5X ," Di sp 1 acemen t  (Tons)  -",F10.2,/ 

+  5 X," St  rut  Length  (Ft.  )  =",F10.2,/ 

+  5X,"S  Wet  Surf  ( F  t  .  *  *  2  )  =",F10.2,/ 

+  5  X  ,  “ L  Center  Buoy.  (Ft.)  - ” , F 1 0 .  2 )  '  )  D I SPT , XST  ,  AREAS , CB 

WRITE! 1 ,'(/ ,5X,  "#  of  Pos.  Wave  Nos.  =",110,/ 

+  5  X  ,  "  2  Hulls?  (  Y  --  1  /  N  =  0  )  :  "  ,  1 1  0  ,  /  )  *  )  N_PO  I  NTS  ,  NSYM 


19 


*****  Calculate  Free  Surface  Elevation. 

** 

CALL  KELVI N_WAVE 
**  , 

*****  Store  Data  File  for  Plotting  Routine. 

*  * 

CALL  WRITE_FILE( FILE2 ,  'KELVIN  WAVE' ) 

** 

NXPL  =  NX_PO I NT 
NZPL  =  NZ_POINT 

WRITE!  1  ,'(/,  5X,  "Data  Written  to  File  :  '',A20)')  FILE2 

WRITE! 1 ,  ' ( / , 5X, “NX_POINT  = " , 1 4 , /5X , " NZ_POI NT  =",I4)')  NXPL, NZPL 

GO  TO  10 

** 

*****  sea  state  Calculations. 

** 

200  CONTINUE 
** 

I ERROR  =  0 
CALL  SEA2D! I  ERROR) 

I F ( I ERROR . NE . 0 )  GO  TO  10 

CALL  WRITE_FILE( FILE2 ,' AMBIENT  SEA  WAVE') 

*  * 

GO  TO  10 

*  * 

*****  Linear  Superposition  of  the  Sea  State  and  Kelvin  Wave. 

** 

300  CONTINUE 

CALL  KELVI N_AND_SEA 
GO  TO  10 

** 

S00  CONTINUE 

** 

STOP 

END 
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CHARACTER*!  ANS 


WRITE ( 1 , 1 

(  5/  ) 

'  ) 

WRITE ( 1  ,  ' 

(  10X, 

"  *  *  J 

******** 

****  MENU  ********** 

****"  )  ' 

WRITE! 1  ,  ' 

(  10X, 

"  ** 

*  *  "  ]  ' 

WRITE! 1  ,  ' 

(  10X, 

M  *  * 

Select 

Item  by  Number 

*  *  "  J  ' 

WRITE! 1  , 

(  10X, 

“  *  * 

*  *  “  )  ‘ 

WRITE! 1  ,  1 

(  10X, 

"  *  * 

(  1  ) 

Kelvin  Wave  Calculations 

*  *  "  )  ' 

WRITE! 1  ,  ' 

(  10X, 

'*  *  * 

(  2  ) 

Sea  State  Calculations 

**”)• 

WRITE! 1 ,  ' 

(  10X, 

”  *  * 

(  3  ) 

Superposition  of  (1)  &  (at 

*  *  "  )  ' 

WRITE! 1  ,  ' 

(  10X, 

"  *  * 

(  4  1 

End  of  Calculations 

**-)• 

WRITE! 1  ,  ' 

(  10X, 

"  *  * 

*  *  “  )  ' 

WRITE! 1  ,  ' 

(  10X, 

"  *  * 

Spec l f  y 

Menu  1  tern  $  “  )  '  ) 

READ! 1 , *  ) 

IDMENU 

IF! IDMENU 

.  LT.  1 

.OR 

IDMENU . 

GT . 4 )  GO  TO  10 

IF  ( I DMENU . EQ . 2  1  THEN 
WRITE ( 1 , ' ( ///10X, 

"**★******************************1 


*  *  *  *  *  *  1 


/10X, 

*  *  *  * 

/10X, 

'  *  *  *  * 

/10X, 

*  *  *  * 

/10X, 

'  *  *  *  * 

/10X, 

‘  *  *  *  * 

/10X  , 

1  *  *  *  * 

/10X, 

'  *  *  *  * 

/10X, 

'  *  *  *  * 

/10X  , 

1  *  *  *  * 

/10X, 

1  *  *  *  * 

/10X, 

'  *  *  *  * 

/10X, 

,  ) 

NOTICE  : 

This  task  requires  that  a  Kelvin 
spectral  file  has  been  read  to 
determine  the  parameters  needed 
for  the  surface  grid  ca 1 cu 1  a 1 1 ons . 
If  a  file  is  not  available  you  will 
be  asked  to  specify  the  dimensions 
of  the  physical  space  tor  the  s^a 
state  calculations. 


****** 
*  *  * 


END  IF 


IF  ( IDMENU . EQ. 3 )  THEN 
WRITE ( 1 , ' ( ///10X, 


******************************************* 

*  *  *  * 

/10X,  ' 

*  *  *  * 

NOTICE  : 

*  *  * 

/10X, ' 

*  *  *  * 

This  task  requires  that  the  Kelvin 

*  *  * 

/I  OX  ,  ' 

*  ★  *  * 

wave  and  the  sea  state  elevations 

*  *  * 

/tox  ,  ' 

*  *  *  * 

data  files  have  been  created 

*  *  * 

/10X,  ' 

*  ★  *  * 

*  *  * 

/I  OX  ,  ' 

*  *  *  * 

Enter  (Y)  if  the  files  exist  and 

*  *  * 

/10X,  ' 

*  *  *  * 

continue  with  the  calculations. 

*  *  * 

/10X, ' 

*  *  *  * 

Or,  enter  IN)  if  the  tiles  have  not 

*  *  * 

/10X  ,  ■ 

*  *  *  * 

been  created  and  return  to  the 

*  *  * 

/10X  ,  ' 

*  *  *  * 

program  menu  to  create  the  tiles. 

*  *  * 

/10X,  ' 

*  *  *  * 

★  *  * 

/  1  OX  ,  ' 

*  *  *  * 

ENTER  Y/N  :$"!') 

READ! 1 , • ( All 

'  )  “.NS 

IF  ( ( ANS . EQ . 

*  n  ' 

OR.! ANS . EQ . ' N ' )1  GOTO  1 0 

END  IF 


RETURN 


SUBROUTINE  SURFACE  _GR I D l  INDEX ) 


************************************************************** 


«*.  *V *\  * 

v.v.N 


csw 


OMMON  /  BLOCK  1  /  X_POINT(  101  )  ,  2_POINT(  101),  Y_PO  I  NIT  (101, 101) 

OMMONI/ B LOCK 4  /  NX_POINT  ,  N2_POI NT , N Y_POI NT  , DX , DZ  ,GRID_SIZE 
OMMON/ BLOCKS /  VEE , AK , YSYM , DI SPT , XST , ZST , AREAS , CB , N_POI NTS , NSYM 

Calculate  Free  surface  Grid  - 

Default  values  for  NX  &  NZ  points  are  61  &  41. 

The  starting  value  for  the  X-axis  (XST)  is  typically 
equal  to  a  ship  length  and  the  axes  system  has  to  correspond 
to  the  axes  used  to  caclulate  the  Kelvin  spectral  functions. 
Currently  the  axes  origin  is  at  the  forward  perpendicular 
with  the  X-axis  along  the  ship  length  taken  as  negative. 


1  =  3. 

RITEt  1 

'  ( 


EAD ( 1 , 
RITEt  1 
EAD ( 1 , 
RITEt  1 
'  ( 


141593 

> 

"  a*******************************************************1 
•'****  Specify  the  Number  of  Points  tor  the  Surface  ****' 
■■****  Grid  in  both  X-  and  2-  directions.  ****' 

"****  Default  Values  are  NX_POINT=61,  NZ_POINT=41  »***' 
"  ****  NX-POINT  =_.")'  ) 

‘  (  13)')  NX_PO I  NT 

, ' ( "***»  NZ- POINT  -_")') 

'  (  13)')  NZ_PO I  NT 


‘★A*********************! 


I***********************1 


'****  Specify  the  starting  value  in  the  X  direction 
'*»**  for  the  rectanqular  grid  in  (ft).  The  default 
****  value  is  equal  to  the  Ship  Length.  Type  /  for 
'****  the  default  value.  XMIN  -  ")') 


********** 
e  in  the  X 
d  in  (ft). 


********** 


he  Rectangular  Gri 
t  .  I  -  ZMAX  .  -*  ZMAX 


DX  =  XMAX/FLOAT ( NX_POI NT- 1 ) 

ZST  =  ELMAX 

DZ  =  2 . 0 * ZST/ FLOAT ( NZ_PO I NT- 1  ) 

GRI D_SI Z  E  -  2 . 0  *  ZST* ( XMAX  -  XST ) 

DO  N  =  1 , NX_POI NT 

X_POINT ( N )  =  XST  -  FLOAT ( N-l ) *DX 
END  DO 

DO  N  «  1 , NZ_PO I  NT 

Z_POINT(N)  =  -  ZST  +  FLOAT ( N- 1 ) *DZ 
END  DO 


*  * 


RETURN 

END 


SUBROUTINE  KELVIN_WAVE 


* 

****************************************************************** 

*  * 

COMMON / BLOCK 1/  X_POINT( 101), Z_POINT( 101 ) , Y_POlNT (101,101) 

COMMON /BLOCK 2/  CGG( 1000 ) ,CGM( 1000 ) 

COMMON/ BLOCK 3/  U_KELVIN(  100  0 )  , V_KELVIN(  1000 )  ,W_KELVINl 100  0 ) 
COMMON/ BLOCK 4 /  NX_PO I  NT , NZ_PO I  NT , NY_PO 1 NT , DX , DZ , GR I D_S 1 Z E 
COMMON /BLOCK 5/  VEE , AK , YSYM , DI SPT , XST , ZST , AREAS , CB , N_POI NTS , NSYM 
COMPLEX  CGG , CGM , DGZ , DGX , EYE 
COMPLEX  CGGN , CGMN 
DIMENSION  STT ( 101  ) 

EYE  =  ( 0 . 0 , 1  .  0  ) 

DO  I  =  1 , N_PO I  NTS 
FI  =  0.00 

IF(I.GT.l)  FI  =  U_KELV I N ( I )  -  U_KELVI N ( I - 1  ) 

IF(  I  .  LT  .  N_POINTS  )  FI  *  FI  +  (J_KELV  I N  (  I  + 1  )  -  U_KELVIN(I) 

CGG ( I )  =  0 . 5*CGG ( I ) *FI 

CGM ( I  )  =  0 . 5*CGM( I ) *FI 
END  DO 


***  Calculate  Free  surface  Grid. 

INDEX  =  1 

CALL  SURFACE_GRID(  INDEX  ) 

DXN  -  DX 
DZN  =  DZ 
XSTN  =  XST 
ZSTN  =  -ZST  +  YSYM 

***  Calculate  Free  Surface  Heiqhts. 

***  Store  in  Matrix  Y_POINT (  I  ,  J  )  . 

DO  N  -  1 , N_POI NTS 

UKELVIN  =  U_KELVIN ( N ) 

VKELVIN  =  V_KELVIN(N) 

CGM(N)  =  CGM(N)*CEXP(-AK*EYE*(XSTN*VKELVIN  -  ZSTN*UKELVI N ) ) 
CGG(N)  =  CGG (N)*CEXP( -AK*EYE*(XSTN*VKELVIN  +  ZSTN *UKEL V I N ) ) 
END  DO 

DO  I  =  1 , NX_POI NT 

DO  J  =  1,NZ_ POINT 
STT(J)  -  0.0 
END  DO 
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DO  N  =  1 , N_POI NTS 
CGGN  =  CGG(N) 

CGMN  *  CGM(N) 

AS  »  REAL (CGGN  +  CGMN) 
BS  =  AIMAG ( CGGN  -  CGMN) 


UKELVIN 

VKELVIN 


U_KELVIN(N) 
V_KELVIN ( N ) 


DGX  =  CEXP ( -AK*EYE*DXN* VKELV I N ) 
CGG(N)  *  CGG ( N ) *DGX 
CGM(N)  =  CGM ( N ) *DGX 

DG2  =  CEXP ( -AK*EYE*D2N*UKELV I N ) 
SD  =  REAL ( DG2 ) 

TD  =  A I MAG ( DG2  ) 

DO  J  =  1 , N2_PO I NT 

STT(J)  -  STT(J)  +  AS 
AST  =  AS*SD  -  BS*TD 
BS  =  BS*SD  +  AS*TD 
AS  =  AST 
END  DO 

END  DO 

DO  J  =  l,N2_POINT 

Y_POI NT ( I , J )  =  STT ( J ) *AK/VEE 

END  DO 

END  DO 

RETURN 

END 
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c  *  *  * 
* 

SUBROUTINE  READ_F I LE ( F I LE_NAME ) 

* 

**************************************************************** 

*  * 

COMMON/BLOCK1 /  X_POINT( 101 ) ,2_P0INT( 101 ) ,Y_P0INT( 101 , 101 ) 

COMMON /BLOCK 2/  CGG( 1000 ) ,CGM( 1000 ) 

COMMON/ BLOCK 3/  U_KELVIN( 1000 ) , V_KELV I N ( 1 000 ) , W_KELVI N ( 1 000 ) 
COMMON / B L0CK4 /  NX_PO I NT , N2_P0 1 NT , N  Y_P0 1 NT , DX , D2 , GR I D_S 1 2  E 
COMMON/ BLOCKS /  VEE , AK , YSYM , DI SPT , XST , 2ST , AREAS , CB , N_POI NTS , NSYM 
COMPLEX  CGG , CGM 
CHARACTER* 20  FILE_NAME 
CHARACTER* 50  TITLE 


V 

'£±& 


; 


K 


100  WRITE(1, 

+  '  ( 

+  ,/,"****  Specify  the  Name  of  the  Data  File  from  which 

+  ,/,"****  this  Program  will  Read  the  Input  Data 

+  ,/,'•****  File  Name  :$")') 


****  " 
****** 


READ( 1 , ‘ (A20 ) ' )  FILE_NAME 

*  * 

OPEN ( UNIT 

=  52  ,FILE  =  FI LE_NAME , STATUS  = ‘OLD1 

*  * 

READ 

(52) 

VEE, AK, YSYM 

READ 

(  52  ) 

DI SPT , XST , AREAS , CB 

READ 

(52) 

N_POINTS , NSYM 

READ 

(52) 

(U_KELVIN( I ) , I = 1 ,N_POINTS) 

READ 

(52  ) 

( V_KELVIN( I ) , I =1 , N_POINTS ) 

READ 

(52) 

(W_KELVIN( I ) , 1=1 , N_POI NTS ) 

*  * 

DO  I 

=  1 

,  N_POINTS 

*  * 

***** 
***** 
***** 
*  * 


READ  (52)  CGG ( I ) , CGM (  I  )  , FW , FFW 
END  DO 

The  varibles  FW  and  FFW  are  not  used  by  this  program  but  they  are. 
These  variables  should  be  left  in  the  above  statements 
for  file  compatibility. 

CLOSE  ( UN  IT  =  5  2 ) 

RETURN 


200  WRITEt 1, ' (4/,5X, 


/,5X ,“****  ERROR 
/  ,5X, "**** 

/  ,  5X  , "**** 
/,5*,“i 
GO  TO  100 

END 


(************" 

; 

-->  File  error  occured  in  opening  the  ****", 
specified  file.  File  does  not  exist.****", 
Tryagain.  ****“, 

r***************"  )  '  ) 
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SUBROUTINE  WRITE_F I LE ( F I LE_NAME , TITLE ) 


COMMON /BLOCK 1 /  X_POINT( 101), Z_POINT( 101 ) , Y_POINT( 101,101) 

COMMON /BLOCK4/  NX_POI NT , NZ_POINT , NY_POI NT , DX , D2 , GRI D_S I ZE 
CHARACTER*  20  FI LE_NAME 
CHARACTER* (* )  TITLE 

100  WRITE ( 1 ,  '  (">>>>>  " , A 2 0 , "  OUTPUT  FILE  <<<<<")•)  TITLE 
WRITE(1, 

+  '(“A*******************************************************" 

+  ,/,"****  Specify  the  Name  of  the  Data  File  For  Storing  *»**" 

+  ,/,“****  The  Output  Results  for  Plotting  ****" 

+  ,/,"****  File  Name  :$")') 

READ ( 1 , ' ( A20 ) ' )  FI LE_NAME 

OPEN ( UNIT  =  5  3 ,FILE=FI LE_NAME , STATU S= ' NEW' , ERR* 200 ) 


WRITE  (53)  ( X_POINT (I) ,1=1, NX_P0I NT ) 

WRITE  (53)  ( Z_POINT ( J ) ,  J  =  1  , NZ_POINT ) 

DO  I  =  1 , NX_POINT 

WRITE ( 53 )  ( Y_POINT ( I , J ) , J = 1 , NZ_POI NT ) 

END  DO 

CLOSE (UNIT*53 ) 

RETURN 


200  WRITE ( 1 , ‘ (4/,5X, 

'•************ 
./,5X,"****  ERROR 
.  /  ,  5  X  ,  "  *  *  *  * 
./,5X,"**** 

.  /  ,  5  X  , "  ************ 
GO  TO  100 


-->  File  error  occured  in  opening  the  ****''  , 
specified  file.  File  already  exists.****", 
Try  again  .  ****''  > 


t* 


Mi 


CLOSE (UNIT*54  ) 

IF(XMIN.NE.X_POINTd  )  )  GO  TO  1010 
IF(ZMIN.NE. Z_POINT ( 1 ) )  GO  TO  1020 
I F ( XMAX . NE . X_POI NT ( NX_POI NT  ) )  GO  TO  1030 
I F ( ZMAX . NE . Z_POINT ( NZ_POI NT ) )  GO  TO  1040 

CALL  WRITE_FILE(FILE_NAME, ‘ KELVIN/SEA  WAVE') 
GOTO  9999 

CONTINUE 
WRITEd, ' (4/,5X, 


/ ,  5X ,  1  ****  ERROR  — >  File  error  occured  in  opening  the  ****''  , 
/,5X,"****  specified  file.  Try  again.  ****", 

/  ,  5  X ,  " ********************************************************" ) •  ) 

IF ( I FLAG . EQ . 1 )  GOTO  100 
IF ( I  FLAG. E Q.2)  GOTO  200 

WRITEd,  '  (4/.5X, 


/,5X, ****  ERROR  — >  XMIN  Does  not  match  for  the  grid  ****"  , 
/  ,  5  X  ,  "  *  *  *  *  from  the  KELVIN  and  SEA  STATE  files.  ****", 

GO  TO  9999 
WRITEd,  '  (4/,SX, 


/,5X,  ****  ERROR  — >  2MIN  Does  not  match  for  the  grid  ****•'  , 
/,5X, "****  from  the  KELVIN  and  SEA  STATE  files.  ****", 

GO  TO  9999 
WRITE(1, * (4/,5X, 


/,5X,  '****  ERROR  — >  XMAX  Does  not  match  for  the  grid  ****'•  , 
/,5X, ''****  from  the  KELVIN  and  SEA  STATE  files.  ****", 

GO  TO  9999 
WRITEd  (4/.5X, 


‘.y/.s'.V 


% 


m 


V  V 


tk  »r 


********************************************************************* 
*  * 

SUBROUTINE  SEA2D( J ERROR) 

*  * 

********************************************************************* 
**  ** 

*****  2-D  SEA  STATE  =  1-D  SEA  STATE  *  COSINE  SQUARED  SPREADING. 

***** 

*****  l -D  SEA  STATE  CALC  USING  ( 1 ) P I ERSON-MOSKOWITZ , 

*****  ( 2 ) JONSWAP ,  ( 3 ) NEUMANN,  (4)  BRETSCHNEI DER . 

*  * 

COMMON/BLOCK1 /  X_POINT( 101 ) , Z_POINT( 101 ) , Y_POINT( 101 , 101 ) 
COMMON/ BLOCK4 /  NX_POINT , NZ_POI NT , NY_PO I NT , DX , DZ , GRI D_S I ZE 
COMMON/BLOCK6/  WM( 101 ) ,AZD(101 ) ,W2D( 101) ,ENF( 101 ) ,EN( 161 ) 
COMMON/BLOCK7/  NSP , FNXP , ISEED , WVKT , THSTD , TZ , NF , I EDW 
DIMENSION  FR( 202 ) ,FS( 161 ) 

** 

CALL  SEA_MENU 

** 

NI  =  10 
TZ  >  10.0 

WVFT  =  1 . 6878*WVKT 
HS  =>  0 . 0182*WVKT**2 
GR  =  32.174 
PI  -  3.141593 
TWOPI  =  2 . *PI 
OMZ  =  0.0 
GAM  »  0.0 

** 

GOTO  (100,200,300,400),  NSP 

** 

*****  p I ERSON-MOSKOWITZ  SPECTRUM  Equation  9,  Reference  (1) 

** 

100  A  *  0 . 0081*GR**2 

B  «  0 . 74*GR**4/WVFT**4 
IA  =■  5 
IB  =•  4 
GO  TO  500 

** 

*****  JONSWAP  SPECTRUM  Equation  10,  Reference  (1) 

** 

200  A  =■  0 . 0  081  *GR** 2 

B  =  0 . 74*GR**4/WVFT**4 
IA  -  S 
IB  =  4 

OMZ  -  ( GR/WVFT )*(C.592)**(0.25)  !  Peak  frequency 

GAM  “3.3  !  Mean  enhancement  factor 

GO  TO  500 

*  * 

*****  NEUMANN  SPECTRUM  Equation  11,  Reference  (1) 

** 

300  A  -  33. 1 *PI /4 . 

B  *  2 . *GR**2/WVFT** 2 
IA  «  6 
IB  -  2 
GO  TO  500 


Equation  10,  Reference  (1) 


!  Peak  frequency 
!  Mean  enhancement  factor 


Equation  11,  Reference  (1) 
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by  Simpson's  rule. 


0 . 000001 


9  * ( 2 . *  F I  -  1  .  )*TWOPI/FN 

/FN 

/FN 

-B/WMI**IB) 

-B/WMR*  *  I B ) 

-B/WML** IB ) 

20 

O  5  20 

0 . 0162*OMZ**2 

0 . 016  2  *OMZ  *  *  2 

0 . 016  2  *OMZ  *  *  2 

(WMI  -  OMZ )**2/bEMI  ;  ) 

(WMR  -  OMZ ) *  *  2 /DENR )  ) 

(WML  -  OMZ ) **2/DENL  )  ) 


550  WRITE! 1 ,'(">>>>>>>>> >  ERROR  <<<<<<<<<<”!') 

WRITE(1,'("<  Intervals  in  intearal  for  E  exceeded  321  > " ) ' ) 

WRITE (  1  ,  '  (  "  <  In  subroutine  SE/.^D  >")') 

I  ERROR  =  3  21 
RETURN 

560  JI  -  J 

FNF  =  NF 

***  Compute  ASW(IW)  and  WM(IW}  based  on  equal  enerqy  intervals. 

DEN  =  0.92*EN(N2)/FNF 

NF1  =  NF  +  1 

ENP  =  0 . 05*EN ( N2 ) 

JE  =  1 

DO  LE  =  2,N2 

570  I F ( EN ( LE )  -  ENP)  590,580,580 

580  FR(JE)  =  FS(LE-l)  +  ( ( FS ( LE ) -FS { LE- 1 ) ) / ( EN ( LE ) -EN ( LE- 1 ) ) ) * 

1  (ENP  -  EN1LE-1  ) ) 

ENP  »  ENP  +  DEN 
JE  »  JE  +  1 

I F ( JE.GT.NF1 )  GO  TO  600 
GO  TO  570 
590  CONTINUE 
END  DO 

***  Compute  E,E2,E4,EP  based  on  FLL  and  FUL  cutoffs. 

600  IF( IEDW.GE. 1 )  GO  TO  610 
DO  IW  -  1,NF 

WM(IW)  =  PIMFR(IW)  +  FR(  IW  +  1)} 

AZDCIW)  *  SQRT(DEN) 

END  DO 

GO  TO  660 

610  DELW  »  TWOPI * ( FR( NF1 )  -  FR(1))/NF 
WL  =  TWOPI *  FR ( 1 ) 

JE  -  1 
FSP  =  FR( 1 ) 

WP  =  WL 

DO  LE  =  2,N2 

620  I F ( FS ( LE )  -  FSP)  640,640,630 

630  ENF(JE)  =  EN ( LE- 1  )  +  ( ( EN ( LE  ) -EN ( LE- 1 ) ) / ( FS ( LE ) - FS ( LE- 1 )  )  ) * 
1  ( FSP  -  FS ( LE- 1 ) ) 

FSP  *  FSP  +  (FR(NFl)  -  FR(1))/NF 
JE  *  JE  +  1 

IF( JE.GT.NF1 )  GO  TO  650 
GO  TO  620 
640  CONTINUE 
END  DO 
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650  DO  IW  -  1,NF 

WN  *  WP  +  DELW 
WM(IW)  =  0 . 5  * ( WP  +  WN) 

WP  =  WN 

AZD(IW)  =  SQRT ( ENF ( IW+1 )  -  ENF(IW)) 
END  DO 

k 

660  DO  IW  =  1,NF 

AZD(IW)  =  1.41421 36*AZD ( I W ) 

WZD(IW)  «  WM(IW) 

END  DO 

k 

CALL  COSQSP ( NF , I ERROR ) 

k 

RETURN 

END 


k  *  * 


SUBROUTINE  COSQSP ( NF ,  1  ERROR ) 

k  a 

AArAAflrAA^frAAAAAAAAAAAAAAAAAAAAAAAAAAAAAAAAAAAAAAAAAAAAAAAAAAAAAAAAAAA 

A  A 

***  2-D  SEA  STATE  CALC  USING  (2*N/PI)*COS( N* THETA )  *  *  2  . 

COMMON /B LOCK  1/  X_POINT( 101  )  , Z_POINT( 101 ) , Y_POINT( 101,101) 

COMMON /BLOCK 4/  NX_POI NT , NZ_PO I NT , NY_PO I  NT , DX ,  "  .  ,GRID_SIZE 
COMMON/ BLOCK6/  WM(101),AZD(101),W2D(101)  , ENF  l 101  )  ,EN( 161  ) 

COMMON /BLOCK 7/  NSP, FNXP, I SEED , WVKT , THSTD , TZ , NF0 , I EDW 
DIMENSION  TR( 101 ) ,TH( 101 ) 

NI  «  10 

EPS  =  1.0E-8 
GR  =  32.174 
PI  =  3 . 141593 
PIRAD  =  PI/180.0 
TWOPI  =  2 . *PI 
NTH  =  9 

I F ( FNXP . GE . 40 )  NTH  *  1 
IF(NTH.EQ.l)  FNXP  -  1 
NAF  =  NF 

***  Compute  the  Energy  Integral  by  Simpson’s  Rule. 


SUMP 
NINT 
J 
N2 
FN 
TH  (  1  ) 
EN  (  1  ) 
SUM 
DO  I 
FI 
TMI 
TML 
TMR 
FNI 
FNR 
FNL 
SUM 


0 . 0 
NI 
1 

NINT/ 2 
NINT 

-  .  5*PI /FNXP 
0 . 0 
0 . 0 
1  ,N2 
I 

( 2 . *FI  -  1 . 


TMI  =  ( 2 . *FI  -  1.  -  . 5*FN)*P1/(FNXP*FN) 

TML  =  ( 2 . *FI  -  2.  -  .5*FN)*PI/(FNXP*FN) 

TMR  =  ( 2 . *FI  -  . 5*FN  ) *PI/( L  NXP*FN ) 

FNI  =  COS( FNXP*TMI ) **2 
FNR  =  COS ( FNXP* TMR ) *  *  2 
FNL  =  COS  CFNXP*TML)**2 

SUM  -  SUM  +  0 . 666667* ( FNR  +  4.*FNI  +  FNL  ) /FN 
EN ( I + 1 )  =  SUM 

TH( 1+1 )  -  TMR 

END  DO 

ER  =  ABS ( ( SUM  -  SUMP ) /SUM ) 

I F ( ER  -  0.001)  130,110,110 

110  NINT  --  2  *N  I  NT 

IFtNINT.GE. 321 )  GO  TO  120 
SUMP  =  SUM 
J  *  J  +  1 
GO  TO  100 
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120  WRITE (1, ' (“>>>>>>>>>>  ERROR  <<<<<<<<<<")') 

WRITE(1,'(“<  Intervals  in  integral  for  E  exceeded  321  >")') 

WRITE ( 1 ,  '  (  “  <  In  subroutine  COSQSP  >")') 

I ERROR  =  321 
RETURN 

*  * 

130  JI  =  J 

** 

*****  Compute  AZD(IW)  and  TM(IW)  based  on  equal  energy  intervals. 

*  * 

N2  =  N2  +  1 

DEN  =  (1.0  -  EPS ) *EN( N2 ) /NTH 

NTH1  =  NTH  +  1 
ENP  =  0.0 
JE  =  1 

*  * 

DO  LE  =  2,N2 

140  I F ( EN ( LE  )  -  ENP)  160,150,150 

150  TRUE)  =  TH(LE-l)  +  ((TH(LE)  -  TH  (  LE-1  )  )  /  (  EN  (  LE  )  -  EN(LE-l)))* 
1  ( ENP  -  EN( LE-1  )  ) 

ENP  =  ENP  +  DEN 
JE  *  JE  +  1 

IF ( JE . GT . NTH  1 )  GO  TO  170 
GO  TO  140 
160  CONTINUE 
END  DO 

*  * 


170  INDEX  =  2 

CALL  SURFACE_GRID( INDEX) 

DO  IX  =  1 , NX_POINT 
DO  12  -  1 , NZ_POINT 
JT  =  0 

IE2  *  2*1 SEED 
CALL  SSEED(IE2) 

ELEV  =  0.0 

DO  IF  =>  1  ,  NAF 

AK  *  W2D ( IF)**  2 /GR 
AZ  *  AZD(IF) 

DO  IW  =  1 , NTH 
JT  =  JT  +  1 
EPT  =  URAN ( JT ) *TWOPI 

TEF  =  0 . 5* ( TR ( IW ) +TR ( XW+ 1 ) ) +THSTD*PIRAD 
AKX  =  AK*COS ( TEF ) 

AK2  =  AK*SIN ( TEF ) 

COSXY  =  AKX*X_POINT ( IX )  +  AKZ * Z_PO I  NT ( I Z )  +  EPT 
ELEV  *  ELEV  +  AZ *SQRT ( DEN ) *COS ( COSXY ) 

END  DO 
END  DO 

Y_PO I  NT ( I X , I Z )  =  ELEV 

END  DO 
END  DO 

RETURN 

END 


35 


SUBROUTINE  SEA  MENU 


COMMON /BLOCK 7/  NSP, FNXP , 1  SEED , WVKT , THSTD , TZ , NF , I EDW 


WRITEC 1 ,'(/)' ) 

WRITE ( 1 **  Specify  the  Spectrum  Type 

,/,"**  for  the  Sea  State  Calculations 
,/,"**  Input  (1)  for  PI ERSON-MOSKOW I  TZ 

,/,"**  (2)  for  JONSWAP 

,/,"**  (3)  for  NEUMANN 

,/,"**  (4)  for  BRETSCHNE I DER 

,/,"**  Specify  the  Appropriate  Spectrum 


$")•) 


READCl,*)  NSP 


WRITE (  1  ,  '  (  /  )  ' ) 

WRITE (  1  ,  '  (  "  *  *  For  60  Deg.,  Spreading  index  =  1.50  **")') 

WRITE ( 1 , ' ( " **  For  90  Deg.,  Spreading  index  =  1.00  **")') 

WRITE ( 1 , ’  ( “ *  *  For  120  Deg.,  Spreading  index  =  0.75  **")') 

WRITE (1,'(“**  Specify  the  Spreading  Index  (ge  .5)  :  $")') 

READCl,*)  FNXP 

WRITEC 1 ,'(/)') 

WRITE ( 1 ,'("* *  Specify  the  Seed  Number” 

,/,"**  for  Random  Phase  Angle  (ISEED.GE.l)  :  S")') 

READCl,*)  ISEED 

WRITEC 1 ,'(/)' ) 

WRITEC 1 , ' ( " **  Specify  Wind  Speed  [S.S.  3  =  16  Kts]=  $")') 

READCl,*)  WVKT 

IF (WVKT  . GT .  0.0  .AND.  WVKT  . LT .  16.5)  TZ  =  7.5 


IF  (WVKT 


WRITEC 1 ,'(/)') 

WRITE ( 1 ,'(“ **  Specify  the  Wind  Direction  " 

,/,"**  with  respect  to  the  X-axis  in  deg.  :  $")') 

READCl,*)  THSTD 
IFCNSP.LT. 4)  GO  TO  50 


0 . 0 

.AND. 

WVKT 

.  LT  . 

16.5) 

TZ 

-  7.5 

16.5 

.  AND  . 

WVKT 

.  LT  . 

21.5) 

TZ 

=  8.8 

21.5 

.  AND. 

WVKT 

.  LT. 

27.5) 

TZ 

=  9.7 

27 . 5 

.AND. 

WVKT 

.  LT  . 

47.5) 

TZ 

=  12.4 

47 . 5 

.AND. 

WVKT 

.  LT. 

55.5) 

TZ 

=  15.0 

55 . 5 

.  AND  . 

WVKT 

.  LE  . 

63.0) 

TZ 

=  16.4 

63.0 

) 

TZ 

=  20.0 

*  *  *  *  < 
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